Excursion Set Halo Mass Function and Bias in a Stochastic Barrier Model of 

Ellipsoidal Collapse 
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We use the Excursion Set formalism to compute the properties of the halo mass distribution for 
a stochastic barrier model which encapsulates the main features of the ellipsoidal collapse of dark 
matter halos. Non-markovian corrections due to the sharp filtering of the linear density field in 
real space are computed with the path-integral technique introduced by Maggiore & Riotto [20]. 
Here, we provide a detailed derivation of the results presented in [22| and extend the mass function 
analysis to higher redshift. We also derive an analytical expression for the linear halo bias. We find 
the analytically derived mass function to be in remarkable agreement with N-body simulation data 
from Tinker et al. [10] with differences < 5% over the range of mass probed by the simulations. The 
excursion set solution from Monte Carlo generated random walks shows the same level of agreement, 
thus confirming the validity of the path-integral approach for the barrier model considered here. 
Similarly the analysis of the linear halo bias shows deviations no greater than 20%. Overall these 
results indicate that the Excursion Set formalism in combination with a realistic modeling of the 
conditions of halo collapse can provide an accurate description of the halo mass distribution. 



I. INTRODUCTION 

Dark matter (DM) is an essential ingredient of the 
Standard Model of Cosmology. Observations provide 
strong evidence that about 90% of the matter in cos- 
mic structures consists of such an invisible component 
[i|-1J] . Central to the DM paradigm is the idea that initial 
DM density fluctuations grow under gravitational insta- 
bility fostering the collapse of baryonic matter. At late 
time, the gravitational infall becomes highly non-linear 
and DM particles violently relax into virialized objects, 
the halos. These are the building blocks in which cool- 
ing baryonic gas falls in to form the stars and galaxies 
that we observe today. Hence, the study of the halo mass 
distribution is of primary importance in Cosmology and 
an accurate modeling of the halo mass function has be- 
come essential to disclosing the complex mechanisms of 
the cosmic structure formation as well as probing the 
physics of the early universe. For example, recent stud- 
ies have focused on the computation of the mass function 
for non- Gaussian initial conditions which are a direct sig- 
nature of inflationary physics (see e.g. [a-0|)- 

The upcoming generation of galaxy cluster surveys (see 
e.g. [8, 9]) will directly probe the distribution of massive 
halos, thus providing cosmological constraints comple- 
mentary to those inferred from measurements of the Cos- 
mic Microwave Background (CMB) radiation and cosmic 
distance indicators. 

The properties of the DM halo mass distribution have 
been mainly investigated using high-resolution numeri- 
cal N-body simulations. These studies have determined 
the halo mass function to a few percent accuracy level 
and provided insights on its redshift and cosmology de- 
pendence |10l - [l3| . In contrast, the development of purely 
theoretical studies has lagged behind. To date, a robust 
theoretical description of the halo mass function is still 
missing. For instance, we do not have a complete under- 



standing of the relation between the conditions that lead 
to the formation of halos and the form of the mass func- 
tion as obtained from N-body simulations. The numer- 
ical analysis usually limits to providing fitting formulae 
which depends on several ad-hoc parameters. 

The seminal work by Press & Schechter (PS) [14!| is the 
first to attempt a derivation of the halo mass function 
from the statistical properties of the initial dark matter 
density fiuctuation field. The idea is that halos form in 
regions in which the linearly extrapolated density field, 
smoothed on a given scale, lies above a given critical 
threshold of collapse, such as that predicted from the 
spherical collapse model [l^- The PS formula explicitly 
depends on such a threshold, which introduces an ex- 
ponential cut-off in the high-mass end as confirmed by 
N-body simulation analysis. However, the PS computa- 
tion suffers of an inconsistent behavior when the variance 
becomes arbitrarily large. In this asymptotic regime it 
is natural to expect that the mass fraction tends to 1 
(i.e. the higher the amplitude of the density fiuctuations 
the larger is the fraction of mass in halos). However, in 
the PS case one finds 1/2, thus suggesting that half of 
the mass in halos has been somehow miscounted. This is 
the so called 'cloud-in-cloud' problem (see e.g. [1^1 and 
discussion therein). 

The formulation of the Excursion Set theory by Bond 
et al. [l7| has provided the Press-Schecther approach 
with a powerful mathematical formalism in which the 
computation of the mass function is reduced to solving a 
stochastic calculus problem (for an exhaustive review of 
the formalism see [18]). As shown in [17'] the smoothed 
density fluctuation field behaves as a stochastic variable 
performing a random walk as function of the smoothing 
scale. Then, halos are associated to random trajectories 
which first cross a critical density threshold of collapse. It 
is the requirement of first-crossing which provides the so- 
lution to the 'cloud-in-cloud' problem. Nevertheless, the 



computation of the halo mass function for reahstic halo 
mass definitions has remained a challenging task. In fact, 
in the Excursion Set theory the mass of a halo depends 
on the form of the filter function, with the latter deter- 
mining the behavior of the random walks. For realistic 
mass definitions the associated filters cause the random 
walks to be correlated, consequently the mass function 
can be estimated only through numerical Monte Carlo 
simulations [13, [l3- Because of this, a thorough sys- 
tematic comparison of the Excursion Set mass function 
obtained accounting for such correlations against N-body 
simulation data has never being performed. 

A major breakthrough in this direction has been re- 
cently made by Maggiore & Riotto [20| who have intro- 
duced path-integral techniques to perform an analytic 
computation of the halo mass function for generic filters. 
This allows us to consistently compare the Excursion Set 
model predictions with N-body data as well with astro- 
physical measurements of halo abundances. 

Here, we compute the halo mass function and the linear 
halo bias for a barrier model which captures the main 
features of the ellipsoidal collapse of halos over a large 
range of masses. In this paper we also provide a detailed 
derivation of the results presented in [23l ■ 

The paper is organized as follows. In Section |TT] we 
briefly introduce the Excursion Set formalism, in Sec- 
tion |lTT] we discuss the modeling of the non-spherical col- 
lapse of halos and the computation of the mass function. 
In Scction lTVl we present the calculation of the corrections 
due to the filter function. In Section |V] we discuss the re- 
sults of the comparison with N-body simulation data. In 
Section IVll we present the computation of the halo bias. 
Finally, we present our conclusion in Section I VIII 



II. THE HALO MASS FUNCTION AND THE 
EXCURSION SET FORMALISM 

In the Press-Schecther approach halos form from re- 
gions of the smoothed linear density fiuctuation field 
which lie above a given density threshold. In such a case, 
the number of halos in the mass range [M, M + dM] can 
be written as 



dn 
dM 



/M 



p d log a ^ 
M2 d\ogM 



(1) 



where p is the background matter density and cr(M) is 
the root-mean-square fiuctuation of the linear dark mat- 
ter density field smoothed on a scale R{M) (containing 
a mass M) which is given by 

CT^(M) = S{M) = -^ I dk k^P(k)W^\k,R(M)], (2) 
2vr^ J 

where P{k) is the linear matter power spectrum at red- 
shift z = and W{k, R) is the Fourier transform of the 
smoothing (filter) function in real space. The function 



f{a) in Eq. ([IJ is usually dubbed as 'multiplicity func- 
tion' and encodes all the effects responsible for the for- 
mation of halos. It is given by /(cr) — 2cr^J^(cr^) where 
J^(5) — dF/dS and F{S) gives the fraction of mass el- 
ements in halos with mass > M . The goal of the Ex- 
cursion Set computation is to evaluate J^{S) and infer 
/((t). Hereafter, we will refer to /(a) simply as the mass 
function. 

Let us briefly review the Excursion Set formalism. The 
density perturbation is defined as S{x) = [p(x) — p]/p, 
where p{x.) is the local density at the comoving position 
X. Then, the smoothed density fluctuation field on a scale 
R is then given by 



J(x,i?)= d'yW{\x-y\,R)S{y), 



(3) 



where W{x,R) is the filter function in real space. In 
Fourier space, Eq. ^ reads as 



(5(x,i?) = 



(27r) 



d^kWik,R)S{k)e 



-ik : 



(4) 



where (5(k) is the Fourier transform of <5(x). Bond et al. 
[17| showed that by taking the derivative with respect 
to R, at any point in space (S(x, R) obeys a Langevin 
equation of the form: 



dS_ 
dR 



= C{R), 



(5) 



where 



am 
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is a noise term, whose properties depends on the statis- 
tics of the underlying density fleld. For initial Gaussian 
density fluctuations, 

(J(k)) - and {S{k)S{q)) = {2nfSD{k + <i)P{k), 

where Sd is the Dirac-function, thus implying that 
(CiR)) = and 

1 /■ -, dW dW 

Halos corresponds to those random walks which hit for 
the first time an absorbing barrier whose value is specified 
by a critical density threshold of collapse. In the spherical 
collapse model this is usually denoted by Sc- 

From Eq. ^ it is evident that the nature of the random 
walks depends on the form of the filter function. For the 
time being, let us assume the sharp-k filter. 



W{k,R) ^0{l/R-k), 



(8) 



then by substituting in Eq. ([7]) and after some algebric 
manipulation we obtain 



dS_ 
dS 



miS), 



(9) 



with {r]s{S)) = and {m{Si)miS2)) = Sd{Si - S2). 
Hence, for the sharp-k filter ??(5) is a white noise and 
d performs a simple Markov random walk as function of 
the variance S. 

Let us define 11 the probability of a trajectory to have 
value S at time S. Then, the probability distribution 
associated to trajectories obeying Eq. ^, which start 
at ^(0) = 0^ and are absorbed by the spherical collapse 
barrier B = Sc, is given by the Fokker-Planck 
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as 



Id'^U 



2 d6^ 



(10) 



with initial condition n((5, 5 = 0) = Sd{S) and absorbing 
boundary n((5c, S) = 0. This solution is 



Il{Sc,S,S) 



V2^ 



fll) 



defined for 6 < 6c- Hence, the fraction of volume oc- 
cupied by halos, i.e. the fraction of trajectories which 
have crossed the barrier {6 > 6c) is given by F{S) = 
1 — jj H((5c,(5, 5')d(5. Evaluating the first-crossing dis- 
tribution dF/dS and substituting in the definition of /(ct) 
we finally obtain 



/(^) 



(12) 



that is the 'Extended Press-Schechter' (EPS) mass func- 
tion. In the next Sections we will discuss how a different 
choice of the halo collapse model and the smoothing func- 
tion alter this standard result. 



III. HALO COLLAPSE AND THE DIFFUSIVE 
DRIFTING BARRIER 

The spherical collapse model '1&\ provides a complete 
description of the non-linear evolution of a spherically 
symmetric top-hat density perturbation embedded in 
a Friedman-Lemaitre-Robertson- Walker background. A 
characteristic feature of this model is that the collapse 
does not depend on the initial size of the region, i.e. on 
the initial radius, but only on the amplitude of the initial 
top-hat overdensity. Since the enclosed mass depends on 
the initial radius, then halos will be indistinctly associ- 
ated to regions of the linear density field which lies above 
the same density threshold. 

The work of Doroshkevich J2S] has shown that initial 
Gaussian density perturbations are highly non-spherical 
and approximately triaxial. Hence, the collapse of a ho- 
mogeneous ellipsoid (see e.g. 24]) should provide a far 
better description of the conditions of halo formation. In 
such a case the collapse depends on the initial overdensity 



^ At very large scales S — > and density perturbation (5 — > 0. 



as well as the shear field and thus on the initial size of the 
ellipsoidal region, i.e. on the enclosed mass. Moreover, 
as also shown in [2^, because of the random nature of 
the density perturbation field, the parameters character- 
izing the ellipsoid are random variables themselves with 
probability distributions which depend on the statistics 
of the underlying density fluctuation field. Consequently, 
the main feature of the non-spherical collapse of halos is 
that the critical density threshold is a stochastic variable 
itself. In the language of the Excursion Set this trans- 
lates into a 'fuzzy' barrier performing a random walk 
whose properties depend on the specificity of the non- 
spherical collapse model considered (see e.g. |25l - [27| ). 
For example, Sheth et al. [2^ have shown that in the 
ellipsoidal collapse model the barrier on average evolves 
as (B) = 6c[l + Pis/ 51)^], with /3 = 0.47 and 7 = 0.615. 
This trend reflects the fact that in the small mass range 
(large S) there is a significant shear that opposes the 
gravitational infall, hence a higher density threshold is 
needed for the collapse to occur. Instead, for large masses 
the shear field is less significant, halos are more isolated 
and the collapse is approximately spherical. 

Moving barrier models have been considered in a va- 
riety of case studies. Sheth [29| has introduced a lin- 
ear barrier model to study the halo distribution in Eu- 
lerian space. In [30| barriers with different functional 
dependence on 5* where used to model the formation of 
halos, filaments and sheets, while in [31| a moving bar- 
rier has been introduce to infer the size distribution of 
H II bubbles during reionization. These studies have 
used Monte-Carlo simulations to infer the mass distri- 
bution. Alternatively, Zhang & Hui ^2| have introduced 
an integral-equation method which allows for a recur- 
sive computation of the mass function for generic moving 
barriers without requiring the use of Monte-Carlo simu- 
lations. Nevertheless, all these studies have focused on 
rigid moving barriers rather than stochastic ones. A real- 
istic description of the non-spherical collapse conditions 
of DM halos must necessarely account for this character- 
istic feature that, as we will see here, modifies the EPS 
result in a very distinctive way. Furthermore, the analy- 
sis of numerical N-body simulations has indeed confirmed 
the stochastic barrier hypothesis. For instance, the work 
by Robertson et al. 33;] has clearly shown that the linear 
density fluctuation associated with halos detected in the 
numerical simulations at a given redshift has not a unique 
constant value, but is randomly distributed as function 
of the halo mass, following an approximately log-normal 
distribution with a nearly linearly drifting average. 

Maggiore & Riotto [21| have implemented these fea- 
tures in the Excursion Set theory by assuming a Gaussian 
diffusive barrier with {B) = 5c and {{B-{B)f) = SDb, 
where Db is a constant diffusion coefficient which they 
set to Db ~ 0.3 to reproduce the results of [3j]. Here, we 
improve the modeling of the diffusive stochastic barrier 
by assuming a linearly drifting average, {B) = B{S) = 
6c + PS (see e.g. [23]), which approximates the average 
drift predicted by the ellipsoidal collapse model in [23]. 





FIG. 1: Schematic representation of random walk trajectories 
for 5 (solid line) and the barrier B (dot line) as function of 
the variance. The variable Y = B — 5 measure the distance 
between the two trajectories and therefore performs a random 
walk as well. The trajectories starts at 5(0) — and -B(O) = 
5c, i.e. y(0) = 5c. 



Under these assumptions the barrier obeys the following 
Langevin equation 



dB 

'dS 



P + Vb{S), 



(13) 



where the noise r]B{S) is characterized by {r]B{S)) = 
and {r]B{S)r]B{S')) = DbSd{S - S'). Without loss of 
generality we can assume that 77s (S*) and r]s{S) are un- 
correlated and the stochastic evolution of the system is 
described by Eqs. (I^- ((T^ with absorbing boundary at 
(5(5) — B{S). The system can be reduced to a one- 
dimensional random walk by introducing the variable 
Y = B — 5. In such a case we have 



dY 



(14) 



with white noise ?7(S') — ris{S)+riB{S) such that {'r]{S)) — 
and {^{S)^{S')) = (1 + Db)^{.S - S'). 

In Fig. [1] we plot a schematic representation of the 
random walks. The variable Y measures the distance 
of the trajectory of the density field to the barrier as 
function of the variance. The system starts at 5 = with 
(5(0) = and B{0) = Sa, thus Y{0) = Fq = ^c Random 
walks which have yet to form halos correspond to F > 
0. The absorbing boundary is located at F = and 
trajectories which have collapsed into halos correspond 
to y < 0. From Eq. p^ we derive the Fokker-Planck 
equation for the probability distribution 110(101^15') of 



FIG. 2: Probability distribution of non-collapsing trajecto- 
ries {¥ > 0) for S = 10'°-3 (solid line), lO"'^'^ (dot line), 
10°-^ (short dash line) and 10°-^ (long dash line). Top left 
panel: Extended Press-Schechter (/3 = 0, Db = 0); top right 
panel: diffusive barrier {j3 = 0, Db = 0.3); bottom left panel: 
linearly drifting barrier (/3 — 0.3, Db = 0); bottom right 
panel: diffusive barrier with linearly drifting average (/3 = 0.3, 
Db = 0.3). 



the random walks with F > (see e.g. |34l|): 



dS 



--/? 



dUp 1 + Db d^Up 
dY ^ 2 9^2 



(15) 



where the subscript '0' refers to the fact that the random 
walks which we are considering are Markovian. We solve 
Eq. (jlSp with initial condition Yo = '5c and impose the ab- 
sorbing boundary condition at F = 0, i.e. no(0, 5) = 0. 
By rescaling the variable Y ^ Y = F/vT+~/7b, a fac- 
torisable solution can be found in the form no(F, S) = 
f/(F,5)exp[c(F - c5/2)], where c = /3/y/TTD^ and 
U{Y, S) satisfies a Gaussian diffusion equation which can 
be solved using the image method [35] or by Fourier 
transform. We find 



no(Fo,F,5) 



-(Y-Yo-l3§) r 



y/2TrS{l + Db) 



2S(1 + D_b) 



(y+yp) 

g 2S(1 + Dg) 



(16) 

Here, it is worth remarking that a general analytic so- 
lution to the Fokker-Planck equation with biased diffu- 
sion and absorbing boundary condition does not exists for 
drift terms which are non-linear in 5. As it will be evi- 
dent in the next Section, having an analytical expression 
for the probability distribution of the Markovian ran- 
dom walk greatly simplify the computation of the non- 
Markovian corrections induced by a realistic filtering of 



the linear density field. For this very reason we have 
opted to assume a barrier model with linearly drifting 
avera ge, rather than the ellipsoidal collapse prediction 
from |28| . 

In Fig. [1 we plot no(ro, Y, S) as function oiY at S ^ 
10-°-3 (solid line), 10-°i (dot line), 10°i (short dash 
line) and lO^"^ (long dash line) for different value oi /3 and 
Db such as to give us a qualitative understanding of the 
barrier model parameter dependence. The standard EPS 
result corresponds to /3 = and Db = (top left panel) . 
The case of a diffusive barrier with (3 = and Db = 0.3 
(as in [2l|) is shown in the top right panel, while the case 
of the linearly drifting average barrier with /3 = 0.3 and 
Db = is shown in the bottom left^ . Finally, the case of 
the diffusive barrier with linearly drifting average with 
/3 = 0.3 and Db = 0.3 is plotted in the bottom right 
panel. 

We may notice that overall amplitude of no(lo, Y, S) 
is a decreasing function of S with an increasing skewness 
toward larger values of Y. Since the total number of tra- 
jectories is conserved, this implies that the probability of 
trajectories that do not cross the barrier between Y and 
Y + dY decreases as function of S, while that of those 
which first-cross it increases. This is consistent with that 
fact that in the bottom-up scenario small mass halos are 
more likely to form than large ones. However, at a finer 
level the trend is barrier model dependent. For instance 
in the case of the diffusive barrier (top right panel) we 
can see that the amplitude of no(yo, Y, S) is smaller than 
the standard EPS result (top left panel). Thus, indicat- 
ing that the number of crossing trajectories is higher. In 
contrast, for the non-diffusive barrier with drifting av- 
erage (bottom left panel) we have that Ilo{Yo,Y,S) is 
larger than the EPS prediction. In addition, the peak of 
the probability distribution rapidly shifts towards larger 
values of Y as function of S as opposed to the EPS case. 
Finally, for the diffusing barrier with linearly drifting av- 
erage (bottom right panel) we may notice that the com- 
bined effect of diffusion and drift is to reduce the overall 
amplitude of 110(10,^, "S*) more effectively in the large 
mass range than in the low mass end. 

The first-crossing distribution gives by definition the 
probability -FqIS) of a random walk to cross the barrier 
between S and S + dS, thus we have: 



MS) = - 



-!sl '^^^oiYo^Y^S) 


(3no{Yo,Y,S) 


°° 1 + Db SHo 
2 dY 


t)„ (Sc+fisf 



S^/^^2n{l + DB) 



(17) 



•^ We find that for P = 0.3 the linear drifting average barrier 
B{S) = Sc + PS approximates to better than 10% the pre- 



fer Db = this coincides with the non-diffusing linear 
drifting barrier solution found in [18l|. Then, the Marko- 
vian mass function reads as: 



Mt) = 



a^/l + DB 



2 __(£c+££ir_ 



(18) 



diction of the eUipsoidal collapse model [2S 
-0.6 < log(l/o-) < 0.4. 



over the range 



The competing effect of diffusion and average drift, 
which we have inferred from the qualitative analysis of 
no(lo, y, S), can be seen directly in the form of Eq. p^ . 
The non-vanishing diffusion coefficient has the effect of 
reducing the amplitude of the mass function cut-off, thus 
shifting it toward smaller values of cr. In contrast, the 
drift of the barrier tends to increase the value of the 
threshold as function of a. 

In Fig. [3] we plot fo{a) for the same values of (3 and Db 
shown in Fig. [51 The qualitative trend confirms that in- 
ferred from the analysis of no(loi Y, S). Firstly, we notice 
that with respect to the standard EPS case (solid line), 
the diffusing barrier model (dot line) shows a cut-off at 
a lower value of a. This is because in the presence of dif- 
fusion the condition of collapse has a scatter around the 
spherical collapse threshold which favors the trajectory 
crossing. Consequently, the mass function is larger than 
the EPS result over a larger range of masses. In contrast, 
the linear drifting barrier (short dash line) gives a mass 
function which is suppressed at all masses. This is con- 
sistent with the fact that the barrier is on average higher 
than the spherical collapse threshold, thus is more diffi- 
cult for trajectories to cross the barrier. In the case of a 
diffusive barrier with linearly drifting average (long dash 
line) the combination of the two effects causes the mass 
function to be tilted with respect to the EPS prediction. 
In the next Section we will discuss the modifications to 
Eq. (IT51) due to a different choice of the smoothing func- 
tion. 



IV. NON-MARKOVIAN CORRECTIONS AND 
PATH-INTEGRAL APPROACH 

A. Halo Mass Definition and Filter Function 

The filtering of the linear density fiuctuation field spec- 
ifies the relation between the smoothing scale R and the 
mass M. The volume selected by the window function 
W{x,R) is V{R) = J d^xW{x,R), hence the enclosed 
mass is given by M{R) = pV{R). However, this relation 
is uniquely specified only for a sharp-x filter, W(r^ R) = 
e{r - R), for which M{R) = 4/3npR^. For generic fil- 
ters, the mass definition is ambiguous, since it is defined 
up to a normalization constant. More importantly, in 
the case of the sharp-k filter, W(k,R) = 6(l/R — k), 
the mass remains undefined (see discussion in [20|)- On 
top of this, we should consider the fact that the mass 
definition of N-body halos depends on the halo detection 
algorithm. Thus, for a consistent model comparison with 
numerical simulations, the filtering of the linear density 
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FIG. 3: Mass function for values of /3 and Db as in Fig. (2] 



field should be chosen consistently with the mass defi- 
nition of the halo detection algorithm used to measure 
the N-body mass function. As an example, the Spherical 
Overdensity (SOD) halo finder detects halos as spherical 
regions of radius Ra enclosing a density pA = Ap, where 
A is the overdensity parameter usually fixed to A = 200 
(which is roughly equal to spherical collapse prediction 
of the virial overdensity at z = in LCDM models) . In 
such a case, the halo mass is Ma — 4/37rpAi?^, that is 
equivalent to definition of the sharp-x filter. However, 
random walks generated by this smoothing function are 
no longer Markovian. In fact, the Fourier transform of 
the sharp-x filter is 



Wik,R) 



(kR)- 



:[sm{kR) - (kR) cos(fci?)]. (19) 



and in such a case it is easy to see from Eq (O that S{S) 
is no longer subject to a simple white noise. 

Maggiore & Riotto [20| have shown that the correla- 
tions induced by the filter function can be treated as 
perturbations about the Markovian solution. Let us con- 
sider a Gaussian random walk, the statistical properties 
are entirely specified by the 1 and 2-point connected cor- 
relators. For a Gaussian field these are {6{R))c = and 



{S{R)S{R')), - 



k^dk 
'2^ 



P{k)W{k,R)W{k,R') (20) 



respectively, with P{k) = Ak'^^T^ik) and T{k) is the 
matter transfer function. In the case of a sharp-k fil- 
ter Eq. (Uni) simply reduces to {S{R[S]WR[S']))c = 
min(S', S'). However, Maggiore & Riotto 'Q3, have shown 
that for the sharp-x filter the 2-point connect correla- 
tor can be written as {S{R[S])6{R[S']))c = min(5',S") -I- 



FIG. 4: Correlation function A{S,S') for a sharp-x filter 
(solid line) obtained by numerically integrating Eq. (|20p for 
5"(R=5 Mpc/h) = 1.228 in the range < S < S' , and the 
fitting function kS/S'{S' - S) with k = 0.475 (dash line). 



A{S,S'), where A{S,S') is a symmetric function which 
vanishes at S" = 5" and is well approximated by 
A{S, S') « kS/S'{S' - S) and k is a fitting coefficient. 
Following their computation, we estimate A(S', S') in 
a ACDM cosmology with model parameters set to the 
WMAP-5yr best fit values: Um = 1-f^A = 0.28, h = 0.7, 
CTg = 0.8 and Ug — 0.96. We cornpute the transfer func- 
tion using the CMBFAST code [Hi]. In Fig. g]) we plot 
A{S, S') for S"(R=5 Mpc/h) = 1.228 (solid line) against 
K S/S'{S'-S) for the best fit value k = 0.475. We can see 
that difference between the numerically evaluated corre- 
lation function and the best fitting formula is well within 
a few per-cent level. Note that k — 0.475 slightly differs 
from the value found in [2^. The discrepancy is due to 
the difference between the numerically computed CMB- 
FAST transfer function and the approximate fitting for- 
mula by Sugiyama [33] which has been used in [20] . The 
coefficient k depends on the assumed cosmological model 
and docs not evolve in redshift. As pointed out in 20j, 
a weakly linear dependence of k on the smoothing scale 
R can improve the fit to the numerically computed cor- 
relation function A{S, S"). Nevertheless, such dependece 
can be neglected to first approximation. 

In terms of the variable Y the non- vanishing connected 
correlators read as {Y{S))c = B{S) = Sc + I3S and 
{Y{S)Y{S'))c = (1 + i:'s)min(S',S") + A{S,S') respec- 
tively. Having assumed the barrier to perform Gaussian 
random walks implies that the non-Markovian part of the 
2-point correlator identically vanishes, i.e. Ab{S,S') = 
0. One may wonder whether a term like Ab{S,S') due 



to the filter function siiould also be included. However, 
there is no reason as to why the barrier has to have the 
same filtering of the linear density field, since the two 
smoothing procedures have very different physical mean- 
ings. The latter is related to the halo mass definition, 
while the former specifies the correlation between the 
condition of halo collapse at different scales. The sim- 
plest approximation is to assume that the collapse at a 
scale S is independent of that at S' , which is equivalent 
to having Ab{S,S') = 0, i.e. Gaussian random walks. 
This suggests that a variety of non-linear gravitational 
effects which induce scale-dependent correlations of the 
halo collapse condition can be implemented in such a for- 
malism through the barrier p-point connected correlators 
with p > 2 or using a non-trival smoothing procedure of 
the barrier random walks. 



B. Path-Integral Method 

Hereafter, we will follow the derivation of [201, ^^i^ ex- 
tend the computation of the non-Markovian corrections 
the case of the diffusive barrier model with linearly drift- 
ing average. 

Let us consider the random walk of the variable Y 
over the time interval [0, S] discretized in steps AS — e, 
such that Sk = ke with k — l,..,n. The probability 



distribution of trajectories that start at Yq, end in Yn at 
time Sn and which never cross the barrier is given by 



^e{Yo,Yn, Sn 

where 



dYi 



dYn_iWiYo,..,Yn,Sn), 
(21) 



W(Yo, ..,y„, Sn) = {SniYiSi) - yi)...fo(y(5„) - Yn)), 

(22) 
is the probability density distribution. Using the Fourier 
transform of the Dirac-function we have 



W{Ya,..,Yn,Sn) 



^:)\e'^U>^or.,^^J-^,x.Y(S^^ 



(23) 

with JVX = JZo^---^- The exponential factor 
within the brackets can be written in terms of the con- 
nected correlators {Y{S))c and {Y{S)Y{S'))c 

(24) 
where we have indicated with Bi = B{Si), Aij = (1 -t- 
DB)min(i, j)e and Ay = nSi/Sj^Sj — Si). Substituting 
Eq. ([H]) in Eq. ([IT]) and given the fact that Ay < 1 we 
can expandi in k, to first order we have 



J 



^<i{Yo,Yn, Sn) 



dYi 



dY„. VX 1 



- z^ij AiAj L^ij 



_^ — lEfc Afc [B}^—Yk\—^TAij\iXjAi 



(25) 



r 



where Aij = (1 + DB)min(i, j). Thus 

Ii,{Yo,Yn,Sn)^in>{yO,yn,Sn)+IV\{Yo,Yn,Sn), (26) 

where H5(io: ^n, •S'n) is the Markovian zero order term 
and Hf (Yoi ^Ti; '5'n) is the first order correction in k. 
The Markovian term obeys the Chapman-Kolmogorov 
equation (see Eq. (jA6p in Appendix |X] for a detailed 
derivation) which is used to explicitely compute the non- 
Markovian term. 

A crucial point concerns the e-dependence of 
Tll^{YQ,YmSn) in the proximity of the barrier. As 
extensively discussed in [20|, the Markovian solution 
Tll{Yo,Yn,Sn) is 0(e) for y„ > and ©(e^/^) at Yn = 0. 
Hence the probability distribution undergoes a transition 
between two different regimes inside a boundary layer of 
finite size. In order to evaluate the form of the proba- 
bility distribution inside this region it is convenient to 
introduce the 'stretch' variable 77 = Y/ ^2e{l + Db) and 
write 

Ill{Yo,Yn,Sn) = C,{Yo,Yn,Sn)u{r^), (27) 

where Cf_{YQ,Yn,Sn) is a smooth solution and u(?/) is a 



function containing the fast variation in e inside the tran- 
sition region. In the continuous limit, lim^_i.oo u^rj) = 1, 
while Cg tends to Eq. pB]) . thus we recover the stan- 
dar d Markovia n solution. Substituting Y = Yn — 
r]yj2e{l + DB) in Eq. ([TBI) and expanding to lowest order 
in e we obtain 



277^0 



^Sl'\l + Db) 



e 2s„(i+db). (28) 



Ge\Y(i, Yn, Sn) 



Hence, substituting Eq. ([28]) in Eq. ([27]) and taking the 
limit 77 — > we obtain 



n^(>^o,o,5„) 



lYo 



Sf/\l + DB) 



-e "'"(i+^B) , (29) 



where 



7= ^ liiR r]u{?]) = -^, 



(30) 



we present the exact derivation of this result in Appendix 
11 
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Similarly we can infer the probability nQ(0,y„,5'„) 
of trajectories starting at the barrier 1^ = and 
ending at Yn > by in troducing the stretch vari- 
able 77 = Fo/-\/2e(l + Db)- Again, substituting Yq — 
r]y/2e{l + Db) in Eq. (IT6l) . expanding to lowest order in 
e and computing Eq. (I27|) in the limit 77 — > we obtain 

S^J\i + Db) 

Finally, the probability of trajectories which start at 
the barrier and end at the barrier can be obtained using 
the dimensional arguments discussed in [201 ■ To lowest 
order in e, we find 

n^(0,0,^„)^e ^ ^. (32) 



C. Non-Markovian Corrections to the Halo Mass 
Function 

We have now all the ingredients to compute the non- 
Markovian correction to first order in k. Using the fact 



that A.e-'-^'t^'-^^'l = -id/dY,e- 
ond term in Eq. p5|) reads as 



iE„AfcfBfc-yt 



'^ ^''l, the sec- 



Ill (1^0 J ^n, ^ri 



IE 



dYi... / dr„-iA,j X 
Jo 



X a,a,Wo(yo,-,i^n,5„), 



(33) 



where Wo{Yq, ..,Yn,Sn) is the Markovian probability 
density of the random walks (see Appendix Ej) . We split 
the sum in Eq. ([33)) as 



1 1 ""^ 

- ^ Aijdidj = -A„„5^ + ^ Aindidn + ^ Aijdidj, 

ij i=l i<j 

(34) 
notice that the first term on the right-hand side of 
Eq. ([M)) vanishes since A„„ — 0, furthermore X]i<7 — 
J27=i J2]=i+i- Then, integrating Eq. ([55)1 by parts and 



using Eq. (IA6|) we obtain 

nu^oi ^ni Sn) = n™i(yo, ^m "^n) + n™i'"(yo, y^, 5„), 

(35) 
where 



n-l 



n™l(^0) i^n, 'S'n) — — 2_^ Aindn 



i=l 



n^(>'o,o,5.)) 



X WoiO,Y„,Sn-S, 



(36) 



is a 'memory'-like term, since it depends on a single sum 
over the past time steps, and 



71 — 2 71 — 1 

n^riYo,Y^,s,,) = j2 E ^^^ 

i— 1 j—i-\-l 



n^(ro,o,5,)x 



(37) 



X n^(o, 0, s, - s,) n^(o, r„, 5„ - s,) 



which represents a 'memory-of-memory' term. For a de- 
tailed derivation of these equations see Appendix O The 
probability distributions in Eq. ([M| and ([57| are given 
by Eq. (P5| . pil) and ([5^ respectively. In order to com- 
pute the sum over the time steps we take the continous 
limit such that 

T,Zi -^ 7 lo dSi and J2t<j -^ ^ lo dS^ J^ dSj. 

We find that the 'memory' term not to contribute to 
the mass function. As a result of the integration over dSi 
we have 



Tl1\Yo,Y,S) 
dY\ 



kYo 



{I + Db^ 

xi-/re^TW(^-^-^t)Erfc 



Yo + Y 



^2S{1 + Db) 



(38) 



and since the first-crossing distribution is given by J^™ ~ 
— -^J^dYIl1^{Yo,Y,S), the subsequent integration of 
Eq. (|58)) over dY vanishes. This is consistent with the 
result of Maggiore & Riotto [2^ . 

The 'memory-of-memory' term cannot be computed 
analytically. Nevertheless, from the ellipsoildal collapse 
model we have /3 < 1. Thus, we can expand the inte- 
grands in powers of (3 and compute the contribution to 
the mass function up to leading order. We find 



£ra—7n{\ 



^5c jia 

(J \ IT 



2 r'2a2 



, (39) 



where a = 1/(1 + Db), k = na and r(0, z) is the incom- 
plete Gamma function. Not surprisingly this expression 
coincides with the K-correction for the diffusive barrier 
obtained in [21|. The first order in (3 reads as 



/i"(^) - -/3«5c [/i'7="('^) + ^Erfc (^^yi 



while the second order is given by 



(40) 






f3 aSciil aSc Erfc 



a 
CT V 2 
1 aS^\ 3 a52 



Sc 



a^ 



4 ct2 



F 0, 



2a2 



2 \ 1 



(41) 



Higher order corrections can be computed semi- 
analytically, since these contain integrals that cannot be 
written in terms of basic functions. Simple numerical 
integration routines are sufficient to evaluate these inte- 
grals. Nevertheless, we find that for /3/(l -(- Db) < 1 
terms 0{> 0^) are negligible. For example in Fig. [5] 
we plot the Markovian mass function /o and the n- 
corrections up to 0(/3'^) and their total contribution for 
/3 = 0.2 and Db = 0.6. We can see that the term 0(/3^) 
is negligible. The largest non-Markovian corrections are 



given by f^aj^{<j) in the intermediate and high mass 
range, and f^~™{a) at low masses. The overall effect of 
these non-Markovian terms is to reduce the overall am- 
plitude of the Markovian solution. 

Notice that the non-Markovian correction f^pJ^{<T) di- 
verges in the very low mass limit (a — >■ oo) due to the 
behavior of the incomplete Gamma function in Eq. p9p . 
As shown in Fig. [SI this term (short dash line) decreases 
to a negative minimum value at log(l/(T) ^ and then 
increases in the low mass range, eventually diverging at 
very small masses, log(l/tT) <^ —0.6 (corresponding to 
M < lO^^h'^M©). This implies that the mass func- 
tion in the case of a diffusive barrier with constant av- 
erage is ill behaved in the very small mass limit. Hence, 
in order to extend its validity one should consider the 
contribution of higher-order corrections in k. However, 
as we have previously discussed, the low mass end of 
the mass function is sensitive to the non-spherical col- 
lapse of halos. In such a case the non-Markovian cor- 
rection due to the average drift of the barrier to leading 
order in (3, /!"I(™(o'), cures the divergent behavior of 

/™^"J2((t). This can be inferred directly from Eq. (HOI , 
where the dependence on the incomplete Gamma func- 
tion is of opposite sign to that in Eq. ([39)) . As shown in 
Fig.[5l /^^"(c) (long dash line) decreases towards neg- 
ative values as /™/3^(o') increases towards positve val- 
ues for small values of log(l/o'). The dependence on the 
Gamma function also appears in /^^™(cr), however this 

correction is subdominant with respect to C(/3). Eventu- 
ally, one can find a small mass range where higher-order 
corrections in k and /3 cannot be neglected. However, 
our analysis shows that in the mass range probed by cur- 
rent simulations, —0.6 < log(l/cr) < 0.4 (corresponding 
10^° < M[h"^MQ] < 10^^) the perturbative expansion 
in K and /? is well behaved provided /3/(l -|- Db) < 1. 



V. EXCURSION SET MASS FUNCTION AND 
N-BODY SIMULATIONS 

The barrier model which we have considered here aims 
to capture the main features of the ellipsoidal collapse of 
dark matter halos. It explicitly depends on /3 and Db 
which parametrize the properties of the collapse thresh- 
old. In principle, these parameters can be determined 
for a given ellipsoidal collapse model. This is because 
the distribution of collapse density values is directly re- 
lated to the probability distribution of the eigenvalues of 
the deformation tensor (see e.g. [2^ [3^). Alternatively, 
one can infer such a distribution by numerically solving 
the ellipsoidal collapse equationsfor randomly generated 
initial conditions (see e.g. 39r;41i]). Then, the values of 
/? and Db can be inferred by best fitting the average and 
the variance of the inferred ellipsoidal collapse density 
distribution. 

Most of the works in the literature have primarely fo- 
cused on determining the average of the ellipsoidal col- 
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FIG. 5: Contributions to the halo mass function /tot (solid 
line) for /? = 0.2 and Db ~ 0.6. The different curves cor- 
respond to the Markovian mass function /o (dot hne) and 
the non-Markovian corrections /i"/3=o (short dash line), /J"g™i) 
(long dash line), /^^fj) (dot-short dash line), f^-^%-, (dot-long 
dash hne). 



lapse threshold [2^, [30|. while no attention has been paid 
to the variance. In l2ll] the authors have provided a 
rough estimate of the variance from the ellipsoidal col- 
lapse barrier numerically determined in [421 ; though they 
did not use such estimate when evaluating the mass func- 
tion. Furthermore, it is very plausible that the distribu- 
tion collapse densities varies with cosmology and redshift. 
Hence, an accurate ellipsoidal collapse model prediction 
of f3 and Db requires a dedicated study which is beyond 
the scope of this paper. 

Here, we test whether the path-integral inferred mass 
function can provide a reasonable description of the nu- 
merical simulation data. In order to perform such a test 
we use the measurements of the halo mass function from 
Tinker et al. [l3| obtained using SOD with A = 200 
on a set of WMAP-lyr and WMAP-3yr cosmological N- 
body simulations. First, we consider the mass function 
measurements at z = 0. For the LCDM models best 
fitting WMAP-lyr and 3yr data the spherical collapse 
model prediction is Sc = 1.673. Using such a value we 
run a likelihood Markov Chain Monte Carlo analysis of 
the excursion set mass function /tot(f) including non- 
Markovian corrections up to 0{/3^) against the data to 
infer the best fit values of /3 and Db ■ The prior parame- 
ter space is log/3 = [—4, 0] and log Db — [—3, 0]. 

We find (3 = 0.057 and Db = 0.294 respectively. The 
numerical simulation data strongly constrain these pa- 
rameters with la errors ajs = 0.001 and aog = 0.001 re- 
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FIG. 6: (Upper panel) Halo mass function at 2; = for the 
Tinker et al. fitting formula with A = 200 (solid blue line), 
diffusing drifting barrier with /3 = 0.057 and Dt = 0.294 
(red dash line) and the corresponding numerical solution from 
Monte Carlo generated random walks with sharp-x filter (red 
squares), Maggiore & Riotto I2l| with Db = 0.235 (green 
dot line). Data points are from [iy|. (Lower panel) Relative 
difference with respect to the Tinker et al. [To] fitting formula. 
The thin black solid lines indicates 5% deviations. 



spectively. We have also verified that these results do not 
change if corrections 0{> (3^) are included in ftot{o')- In 
Fig- m (upper panel) we plot the best fitting mass function 
(red dash line) against the simulation data together with 
the four-parameter fitting formula by Tinker et al. [lO| 
for A = 200 (solid blue line). For comparison we also plot 
the diffusive barrier by Maggiore & Riotto [21] best fit- 
ting the data with Db = 0.235 (green dot line). In Fig. [5] 
(lower panel) we plot the relative differences with respect 
to the Tinker et al. formula. We may notice the remark- 
able agreement of the diffusive drifting barrier with the 
data. Deviations with respect to Tinker et al. (2008) are 
within the ~ 5% level for log(l/o') < and within 7% 
over the range < log(l/o') < 0.3. This is quite impres- 
sive given the fact that our model depends only on two 
physical parameters. As expected the improvement with 
respect to the diffusive barrier [2l[ is due to the drifting 
average which systematically suppresses the formation of 
small mass halos with respect to the massive ones. In 
Fig. [S] we also plot the mass function for the diffusive 
barrier model with /3 — 0.057 and Dj, = 0.294 inferred 
from Monte-Carlo generated random walks with sharp-x 
filtering (red squares). As it can be seen this numerical 
solution well reproduce our mass function formula with 



FIG. 7: Redshift evolution of the mass function residuals with 
respect to the redshift dependent Tinker et al. [lOl fitting 
function (blue solid line) for diffusive barrier with linearly 
drifting average (red dash line). 



the same level of deviations from the Tinker et al. data'^. 

It is worth noticing that the best fit value of /3 is about 
a factor 5 smaller than the ellipsoidal collapse model ex- 
pectation. This could be an artifact of the SOD mass 
function measurements, a consequence of our modeling of 
the barrier diffusion as a Gaussian random walk rather 
than log-normal, or having limited the non-Markovian 
corrections to first order in k. In a future study we will 
perform a more detailed analysis to discriminate between 
these different effects. This will allow us to obtain an 
unbiased physical interpretation of the barrier model pa- 
rameter values. 

The mass function measurements from Tinker et al. 
[10l | extend to z = 0.5, 1.25 and 2.5. The mass function 
can also reproduce these measurements for a given com- 
bination of values of (3 and Db- At these redshifts we 
have (5c = 1.680, 1.685 and 1.686 respectively ^. 

In Fig. [7] we plot the residual of the mass function 
/tot(o') for the barrier model parameters best fitting the 
data at z = 0.5, 1.25 and 2.5 with respect to redshift de- 
pendent fitting formula by Tinker et al. [l3l- We can 



^ Recently, in [43 | the authors claim that the path-integral ap- 
proach does not reproduce the Monte-Carlo inferred mass func- 
tion to no better than 20% and thus conclude that such a for- 
malism is inappropriate to fit halo abundances dn{M)/dM. Our 
numerical evaluation for the diffusive barrier model with linearly 
drifting average shown in Fig. |6] clearly demonstrates this not to 
be the case. 
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see that, apart a systematic deviation at high masses 
(log(l/(T > 0.2)) of order of < 10%, ftoticr) is consistent 
with data to better than 5%. In particular at z = 2.5 
the functional form of the halo mass function seems to 
better reproduce the numerical measurements compared 
to the Tinker et al. [iO| fitting formula for which the 
authors have found residuals > 5% at z = 2.5. In Fig. [5] 
we plot the best fit values of (3 and Db as function of 
z. The value of /3 tends to saturate at z > 0.5, which 
would imply that the non-spherical collapse threshold on 
average has a similar mass (scale) dependence at higher 
redshifts than at present. On the other hand, Db peaks 
at z = 0.5 and then decreases. Since, the diffusion co- 
efficient primarily affect the high-mass end of the mass 
function, this would suggest that the collapse threshold 
of massive halos at higher redshifts is closer to that pre- 
dicted by the spherical collapse. However, because of 
the limited z-sampling as well as a systematic bias of 
the mass function data towards high masses the physical 
interpretation of these trends should be taken carefully. 
We leave a detailed analysis of these dependencies to a 
future numerical study. 



VI. LINEAR HALO BIAS 

Halos are biased tracers of the dark matter density 
perturbations from which they form. Operationally the 
halo bias is defined as the ratio of the 2-point halo spatial 
correlation function to that of the underlying dark mat- 
ter density fluctuation field. In the Excursion Set for- 
malism this can be estimated using the peak-background 



split technique (see e.g. [4^-|47[) The basic idea is to 
evaluate the conditional first-crossing distribution and in- 
fer the relative abundance of halos of a given mass (i.e. 
S) as function of the large scale density fluctuation ^Oi 
^"(51^0, S' — 0). Then, it can be shown that to first order 
in Sq the halo bias is given by 



6h(5) 



1 



1 



dHS\So,0) 



J"(5|0,0) dSo 



(42) 



where 7^(5*10,0) coincides with the unconditional first- 
crossing distirbution with Sq = 0. 

In the framework of the Excursion Set theory the halo 
bias for the sharp-x filter with the non-Markovian correc- 
tions to first order in k has been derived in [4a|. We ex- 
tend their calculation to the diffusive barrier with linearly 
drifting average (see also [49[ for the case with sharp-k 
filter and 50] for a computation of the conditional mass 
function in the case of generic moving barrier models and 
sharp-x filter). The calculations are quite cumbersome 
and since the basic results by Ma et al. [45 applies also 
to our case, we will report only the relevant passages. 

In order to evaluate ^"(51150, 0) let us first compute the 
conditional mass function with conditioning on a generic 
scale S' < S where d{S') < B{S'). We find convenient 
to work with the variable Y = B — 6, even though the 
conditional first-crossing distribution T{S\6',S') differs 
from J^(S'|y', S"). In fact, the latter imposes the condi- 
tion on both variables S and B, while the former does not 
impose any condition on the barrier value. Nevertheless, 
since we are interested in computing the first-crossing 
distribution in the large scale limit, the barrier trajecto- 
ries converge toward a unique constant value, Bq — Sc 
for S" — ?> and thus we recover J^{S\So,0). Nevertheless, 
some care is needed when computing the non-Markovian 
corrections. 

Following 48] the path-integral defition of the condi- 
tional first-crossing distribution is 



/•oo 
■J (i^nl^mi ^m) ^ I din 

Jo 



O-T (l„, bn\irm i>m) 



(43) 



where 



j;^ dYi...dY^_iW{Yo = Sa,...,Y^,S^) 



(44) 



with the probability density developped to first order in 
K reads as 



n 

W{...) ^Wo{...) + - Y. XAdjWoi.. 



(45) 



i-J = i 



and Wo(...) is the probability density distribution of the 
discrete Markovian random walks. Notice that the inte- 
gral in the denominator of Eq. (|44p provides the correct 
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normalization factor to the conditional first-crossing dis- 
tribution. The numerator in Eq. (|44p can be computed 
by splitting the sum in Eq. ()45|) and computing each term 
individually. However, as shown in 48] only few of these 
terms actually contribute to J-'{Sn\Ym, Sm) and we have 
verified this to be the case also for the diffusive barrier 
model with linearly drifting average. In particular we 
have 

7Vf(y„,r„,^„,5„) (46) 



+ -1-1-1 (,^mi ^nj "^n ^m) ~^ 



no(^0 — Sc: ym, Sm) 



where Ilo{Ym,Yn,S„ - Sm) and no(Fo ^ Sc,Ym,Sm) are 



given by the Markovian probability distribution Eq. (|16p . 
The other two terms in Eq. (|46| contains the non- 
Mar kovian corrections to first order in k, these read as 



J 



-1 ''( ~ 1 /"OO />oo 

^l{Ym,Yn,Sn ~ Sm) — 7: 2_^ ^'J / ^Ym+l--- I dYn-ldidjWf){Ym-, ■■■,Yn, Sn — Sri 



and 



■■,3 

n-l 



(47) 



^ A,,n^(K„„ 0, s, - 5„)n^(o, 0, s, - 5.)n^(o, y^,, 5„ - s,), 

2 j'^rn+l 



^~1 /'OO /"CXD 

7Vi^(y,„,y„,5„,5„0- ^ Aj„, / dYi...dYj... dYn^idmWo{5c,:;Ym,Sm.)WQ{Ym,...,Y,,S,^Sm.y, 

0=m+l Jo Jo 

X W^o(Yj, ..., F„, 5„ — Sj) = dm^o{Yo = Sc, Ym, Sm) 2_^ Aj,„nQ(y,„, 0, Sj — S'm)nQ(0, Yn, Sn ~ Sj). 

J— ?n-t-l 



(48) 



r 



As in the case of the non-Markovian corrections to the 
mass function, these terms can be computed in the conti- 
nous limit with the sum sustituted with an integral over 
the variance and the integrands given by Eq. (f29l) . (|3T|) 
and (|32|) respectively. 



The first term in Eq. (1461) gives the Markovian condi- 
tional first-crossing distribution. In the limit Sm — > 
(i.e. Ym -^ 5c — (5o) we find 



^o(^l<5o,0) = 



(^C — 5q) fir ^ (-^0-^0+135)2 



(49) 



Eq. (j47p is equivalent to the 'memory-of-memory' term 
Eq. (|37l) In the limit Sm — > 0. The double integral over 
the variance can be computed analytically by Taylor ex- 
panding in /3. Here, we limit the computations to terms 
up to first order in /3, since its numerically calibrated 
value is « 10~^. We find 






,(£c-io) 
' 23 — 



ir 
2 



27r 
a{5c 



So)' 



2S 



(50) 



and 



J^l^,r>{S\S„,0) 



-/3a((5c-(5o) 
'Sr- 



J-i%o(^l'5o,0) + 



25* 



-Erfc 



25* 



(51) 



On the other hand, the computation of Eq. P5|) re- 
quires some care. This term contains the derivative of 
the Markovian solution, Ilo{6c,Ym, Sm), with respect to 
Ym, hence it explicitely depends on the variation of the 
probability distribution with respect to the distance of 
the barrier Bm to Sm- Nonetheless, we are interested 
in the first-crossing distribution which is conditional in 
Sm only, and for which the result of the integration of 
Eq. pS)) in the limit Sm ~^ should depend only on the 
variation with respect to Sm- However, after computing 
this term we find that for /3 = Db = we do not re- 
cover the non-Markovian correction of spherical collapse 
model. This is a direct consequence of the fact that we 
are working with Y rather than i5. The correct way to 
proceed is to first marginalize over Bm and then take the 
limit Sm -^ 0. However, this is a very cumbersome com- 
putation, instead we have found that the inconsistency 
can be cured by simply taking the derivative in Eq. (j48p 
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with respect to Sm only. In such a case we find 

1 



{Sc - So)' 



(52) 



which for Db = coincides with the resuh of [48| and 



K,pw iS\So.^) = -Pa K -^^ ^ V 2^^ 

{6c - So)' 



X e 



^" — js — 



l-a 



S 



(53) 



Then, summing ah these terms and evaluating Eq. (j42p 
we finally obtain 



bhi:^) = 1 + 



1 



1 



eT-r 



(o,^) 



Sc 1 

/3a 



Ai + |e^r(0,^) 



l_^+|e-#r(0,^) 



1-K 



e-r(o — 




(54) 



where z^ = Sc/(t. In Fig. [S] we plot the halo bias for 
/3 and Db best fittting the mass function data at z = 
(red short dash line) against the best fit formula to 
the bias measurements inferred in [5l| (black solid line). 
As we can see the difference is < 20% over the mass 
range probed by the simulations. This is consistent with 
the findings of Ma et al. [4^. As argued in [51] the 
discrepancy with respect to the halo bias data from N- 
body simulations is related to the peak-background split 
approximation itself. I Hence, the improvement on mass 
function calculation does not give any further insight on 
the linear halo bias. 



VII. CONCLUSIONS 

The Excursion Set formalism provides a powerful 
mathematical framework which allows us to perform a 
theoretical computation of the halo mass function from a 
limited set of initial assumptions. These must involve the 
statistics of the linear density fiuctuation field as well as a 
stochastic barrier model of the halo collapse conditions. 
In addition, such calculation needs to be implemented 
with a path-integral evaluation of the corrections due to 
the filtering of the linear density field associated with a 
realistic mass definition. Such an approach allows for 
a consistent model comparison with N-body simulation 
data. 

Here, we have derived an analytical expression for the 
mass function and linear bias in the case of a diffusive 
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FIG. 9: Top panel: Halo bias for the diffusive drifting bar- 
rier model with parameters best fitting mass function data 
from [10| at jz = (red short dash line) against the best fit 
formula to the halo bias from the same numerical simulation 
sets inferred in [S^] (black solid line). Lower panel: relative 
difference with respect to the halo fit formula from 51]. 



barrier model with linearly drifting average. This model 
well approximates the main features of the ellipsoidal 
collapse. We have found a remarkable agreement with 
N-body simulation data with differences « 5% over a 
large range of masses. Such an agreement is due to the 
competing effects of the barrier average drift at small 
masses and of the diffusion in the high-mass end. This 
has important phenomenological implications especially 
in the study of primordial non-Gaussianity. In fact, sev- 
eral studies have estimated the halo mass function in 
the case of non-Gaussian initial conditions assuming the 
spherical collapse models (see e.g. [3, 0, [521 )• However 
the comparison with non-Gaussian N-body simulations 
has shown large deviations in the low mass range com- 
pared to large masses (see e.g. Fig. 1 in [5J]). In the 
light of our results, it is plausible that such discrepancies 
may be attributed to the non-spherical collapse of halos. 
The inclusion of a simple diffusive barrier model with lin- 
early drifting average for a non-Gaussian linear density 
field could resolve or alleviate the problem. 

Our results suggest a number of directions which war- 
rent further investigation. Firstly, it will be insightful 
to derive the statistical properties of the fuzzy barrier 
for a given ellipsoidal model as function of the variance 
of the linear density field. This will provide theoreti- 
cal predictions for /3 and Db which can be confronted 
with numerically calibrated values for different redshifts 
and cosmologies, and it will allow us to accurately test- 
ing the modeling of the halo collapse condtions. On the 
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other hand, in the upcoming years several observational 
campaigns will probe the halo abundance through galaxy 
cluster surveys. The mass function derived here can be 
used to perform a data analysis of the barrier model pa- 
rameters, thus providing information on the collapse of 
DM halos which has been previously unforseen. 



Finally, from Eq. ((2T|) we obtain the relation 

Jo 

(A6) 

which is the Chapman-Kolmogorov equation for a Gaus- 
sian random walk with linearly drifting average. 
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Appendix A: Chapman-Kolmogorov Equation 

The probability density distribution of the discrete 
Markovian random walk is given by Eq. (1241) which by 
expliciting the variance dependence of the connected 1 
and 2-point correlators reads as 

(Al) 
where Aij = (1 + £'B)min(i, j'). Diagonalizing the 
quadratic form in Eq. (|Aip and solving the resulting 
Gaussian integral we obtain 

[2'!Te{l + Db)]^ 

(A2) 
where Aij ~ min(i,j). One can show through induction 
that (A"i)ij = 2 for i = l,...,n- 1, (A"i)„„ = 1 and 
{A-'^)ii+i = {A-^)i+i^ = -lior i = l,...,n-l, while all 
other elements vanish. Similarly one finds that det A = 1. 
Thus, we can write Eq. (|A2p as 



Wo{Yo,..,Y„,Sn)^ 






= Ve(Ar)VKo(ro,-,i^«-i,5„_i), 
with Ay = Yn — Yn-i and 



(A3) 



M^Y) 



1 _(&Y-0_ef_ 

-.e 2.(i+db) , (A4) 



^27re(l + Db) 

An important consequence of the above relations is 
that the Markovian probability density satisfies the rela- 
tion: 

Wo{Yo,...,Y,^0,...,Yr,,Sn)^Wo(Yo,...,Y,^i,0,S,)x 

xWoiO,Y,+i,...,Yn,Sn-S,) 

(A5) 



Appendix B: Calculation of 7 

We can compute the factor 7 defined by Eq. pO)) sim- 
ply using the properties of the Markovian solution to 
the Fokker-Planck equation Eq. (fT5|) . Without loss of 
generality we can consider the case where the absorbing 
boundary is at F = Yc < rather than Y = 0. In such 
a case the solution to Eq. p5|) with initial condition at 
Y = Yq > Yc reads as 



n 



ga/3(y-yo-/3f) /^ 



V2^ 



2S 



a(1Yc-Y ~YqY 
■ e 2S 



(Bl) 

where a = 1/{1+Db)- Taking the derivative with respect 
to Yr we have 



dn 



2a2Yc-Y(i-Y 



a(2Y^-Yo-Y)^ 



S3/2 



_^al3(Y-Yo-fi^)^-- 5^ 



(B2) 

On the other hand using the path-integral formulation 
we have 



thus 
'dY^^ 



dYi... / dYn^iW{YQ,...,Yn,Sn). (B3) 



E 



dYi. 



/B4) 
then using Eq. (|A5P and taking the continous limit we 
obtain 

|^--Hmi / ds,n,(o,i;,5,)n,(ye,y,5-5,), 

dYc e^O £ Jq 

(B5) 
where the left-hand side is given by Eq. (|B2p . while 
ne(0, Yc, Si) and lieiYc^ Y,S — Si) are the finite-e correc- 
tions to the Markovian solution near the boundary. As 
shown in Section IIVBI these terms can b e evaluated by 
introducing the stretch variable ri^j2t/a and expanding 
the Markovian solution to lowest order in e, this gives 

ne(ro,y.,^)^Vi^^%^e-^^^-^#^, (B6) 



53/2 



and 



ne(re,r,5) = Vi^^^§7^e-^^^^-^^^. (B7) 
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Substituting these solutions on the right-hand-side of then, substituting in Eq. (|B8I) and comparing with 

Eq. (|B2p after simplifications we obtain 7 — 1/ ^/tt. 



(B8) 



Eq. (|B5p we have 
dti 



dYr 



{aj)'{Y,-Yom-Y)e^^^^~''"-^i^x 



dS-, 



ajY-Yg)-' a{Yo-Yay 

g 2(S-Si) g 2^i 



the integral in dSi can be computed analytically since 



£25- 2(S-S,) ^ + 5 1 (A + B)^ 

dS,.^r77^ = V2tt — 7^^-377:76 5^~ 



' sf^^S ~ S^)^/^ 



AB 53/2 



Appendix C: Non-Markovian Corrections 



Using Eq. p4p the non-Markovian correction to first 
(B9) order in k, Eq. ([33|. has two separate contributions 



= -VA™a„ / dY^...d%... f 

Uo Jo 



/ ^ ^inCfn 



dYi...dY,... / dYn_iWo{Yo,..,Y,^0,S,)Wo{Y,^0,..,Yn,Sn-S,) 



n^(ro,o,5,)n;^(o,y„,5„-5,) 



(Cl) 



where we have used Eq. (jASp . and 



-ididj 



/•OO /"OO 

n^i'"(ro,F„,5„)- Va„ / dYi... / dr„_ia, 

"'O Jo 



j<j 



dYi...dYj... / dr„^i9, 
Jo 



Va,j / dyi...df,...dr,... / dr„_i 

^^j Jo Jo 



WoiYQ,..,Yn,Sr 

WoiYo,..,Y]=0,..,Yn,Sn) 
WoiYo, ..Y, = 0, .., Y, = 0, Sj) Wo{Y, = 0, .., y„, 5„ - S,) 






n^(yo, 0, 5,) n^(o, o, s*, - s,) n^(o, r„, 5„ - s,) 
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